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Abstract 

The theory of plastic bending of single crystal beam having two active slip 
systems is proposed. Applying the variational-asymptotic method we reduce 
the energy functional of the beam to the one-dimensional energy functional 
which admits analytical solutions for the symmetric slip systems. The thresh¬ 
old value at the onset of plastic yielding, the dislocation density, as well as 
the moment-curvature curve are found. We solve also a similar problem that 
takes the energy dissipation into account. 
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1. Introduction 


When a single crystal beam is plastically bent, excess dislocations occur 
on certain active slip planes to reduce its energy. As a rule, the number of 
dislocations is huge, so, instead of a difficult task of describing individual 
dislocations, the use of continuum approach to hnd the dislocation densities 
in the bent state turns out more effective. The first attempt of taking into 
account continuously distributed excess dislocations in the plastically bent 


beam was made by Nye (1953) who expressed the curvature of the beam 


caused by dislocations in terms of the dislocation density tensor bearing now 
his name. Read (1957) and Bilby et ah (1958) have extended this result 


to the case when the stress due to excess dislocations does not vanish (see 
also (Ashby, 1970 Evans, 1995 Mura, 1987 Wang et ah, 2003)). Unfortu¬ 


nately, because of the lack of physically sound continuum dislocation theories, 
no quantitative dislocation distribution could be predicted in terms of the 
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bending moment. The first qnantitative prediction of the dislocation density 


of single crystal beam bent in single slip has been made in (Le and Nguyen 


2013). The analyzed model was based on the continuum dislocation theory 


proposed by 

Berdichevsky 

(2006al|b 

) and developed further in (Berdichevsky 

and Le, 2007 

IKaluza and ^ 

4|2011t 

Kochmann and Le, 2008, 2009a|b[ |Le and 


Sembiring 


2000 Gur 


2008a 


m 


b, 2009) (see also the alternative approaches in (Acharya 


2002 Gurtin et ah, 2007)). The information about dislocation 


distribution obtained in ( Le and Nguyenf 2013) has been used to explain 
the formation of small angle tilt boundaries and predict the number of poly¬ 
gons in a polygonized beam after annealing (see also (Le and Nguyen, 2011| 


2012 )). 


The real bending of single crystal beam is often more involved because 
two or more slip systems could be activated during the plastic bending, de¬ 
pending on the orientation of crystal with respect to the beam axis affecting 
the Schmid factors. For instance, if the fee single crystal beam having the 
axis parallel to [Il0]-direction is bent about the [110] direction, then there 
are two symmetric active slip systems oriented at angles ±54.7° with respect 
to the beam axis (see, e.g., (Kysar et ah, 2010|[ )). Simulations of such bending 
problems require more elaborate models because of the interaction between 
edge dislocations of different groups. This paper aims at constructing an 
approximate one-dimensional theory of bending of single crystal beams hav¬ 


ing two active slip systems and, thus, extending the results obtained in (Le 


and Nguyen, 2013) to double slip. We incorporate the interactions of excess 


dislocations in a material model similar to that proposed by |Le and Sem- 


biring (2009). We then use the variational-asymptotic method developed 


in (Berdichevsky, 1983 Le, 1999) to derive the one-dimensional approximate 


theories from the exact three-dimensional problems in two cases: (i) the resis¬ 
tance to dislocation motion is negligibly small so that the displacements and 
plastic slips can be found by minimizing the energy functional, (ii) the rate- 
independent dissipation is taken into account. The obtained one-dimensional 
variational problems admit analytical solutions representing the smooth min- 
imizers for the symmetric active slip systems. It is established that there exist 
energetic as well as dissipative thresholds for the dislocation nucleation which 
depend on the thickness of the beam. This exhibits the typical size effect of 
the gradient theory. Based on this analytical solution the deflection of the 
beam, the dislocation density, and the moment-curvature curve are analyzed 
in terms of the bending moment for different loading/unloading processes. 

The paper is organized as follows. In the next Section the 3-D models of 
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single crystal beam with continuously distributed dislocations bent in double 
slip are formulated. Sections 3 and 4 derive the one-dimensional equations 
and solve them for the case of zero dissipation. Sections 5 and 6 consider the 
same problems but with dissipation. Finally, Section 7 concludes the paper. 

2. 3-D models of crystal beam bent in double slip 



Figure 1: Single crystal beam bent in double slip by a moment M. 


Consider a single crystal beam bent by a moment M. In the undeformed 
state the beam occupies the domain V of the three-dimensional euclidean 
space, with the cartesian coordinates {x,y,z) G (0, L) x (—h/2, h/2) x (0, a) 
where L, h, and a are the length, height, and width of the beam, respectively. 
We will assume that h a ^ L and that the beam is bent under the plane 
strain condition. For dehniteness, let the beam be clamped at x = 0 and 
subjected to a linearly distributed traction —ry having the resultant moment 
M at X = L. Besides, the upper and lower boundaries of the beam aXy = h/2 
and y = —h/2 are traction free. If the applied moment is small, it is natural 
to expect that the beam deforms elastically and the stress distribution is 
linear over the thickness according to the elementary beam theory. However, 
if the applied moment exceeds some critical value, edge dislocations with 
dislocation lines parallel to the z-axis may appear to reduce energy of the 
bent beam. We assume that, at the onset of yielding and during the plastic 
deformations, the crystal admits two active slip systems with slip planes 
inclined at angles (pi and (p 2 to the plane ?/ = 0 as shown in Fig. 

Under the plane strain state condition there are only two non-zero com¬ 
ponents of the displacements that do not depend on Ux{x, y) and Uy{x, y). 
Consequently, the non-zero components of the total strain tensor are 


^XX ^X,X) 


^yy '^y,y> 


^xy ^yx c2^^x,y 4 “ ^y,x)- 
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Throughout the paper the comma before an index is used to denote the 
partial derivative with respect to the corresponding coordinate. Since there 
are two active slip systems, the plastic distortion tensor is given by 


Pij = I3i{x, y)s]m] + /32(x, y)slm], 

with = (cos (pi, sin(pi, 0 ), = (— sin<^i, cos(pi, 0 ), = (cos(p2, sin (p2, 0 ), 

= (—sin 932, cos </92,0) being the unit vectors denoting the slip direction 
and the normal to the slip planes of the corresponding slip system. We call 
Pi^x^y) and (32{,x,y) the plastic slips. Thus, the hnding of plastic distortion 
tensor reduces to determining the functions (3i{x,y) and (32{.x,y). Non-zero 
components of the plastic strain tensor, + /Sjj), read 

= -^/5isin2(y9i - ^/?2sin2v32, ely = ^/3isin2(y9i + ^/32 sin 2(^2, 

COs2(pi + ^^2 cos 2(^2. 

Accordingly, the non-zero components of the elastic strain tensor, efj = Sij — 
^Ix = Ux,x + sin2(y9i -F ^/32sin2(y32, 

^yy = ~ sin2(pi - sin 2(^2, 

=^lx = \ i'^x,y + COs2(pi - COS 2(^2) • 


The distribution of excess dislocations associated with these two active 
slip systems is described by the dislocation density tensor (introduced by 
Nye (1953), Bilby (1955), and Kroner (1955)), aij = ejkil3ii^k, whose non-zero 


components read 


ctxz = /3i,x COS^ (fi + f3i^y cos ifi sin cpi + /32,x cos^ (p 2 + l32,y cos (P 2 sin (p 2 , 
Oiyz = (3i,x cos (pi sin (pi -f sin^ (pi -F /32,x cos ip 2 sin (^2 + /32,y sin^ (p 2 - 


Thus, we see that there are two different groups of excess edge dislocations 
associated with the slip systems 1 and 2, whose resultant Burgers’ vectors 
show in the directions and s^, respectively. These give two scalar disloca¬ 
tion densities 

Pl = ^\Pi^xCOSipi+Pi^ySinipil, P2 = ^\P2,xCOSip2 +P2,ySmip2\, ( 1 ) 
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with b being the magnitude of the Burgers vector. 

Let us propose the free energy density of the single crystal deforming in 
double slip in the form (see (Le and Sembiring, 2009)) 


= In ^ 


1_ PI 

Ps 


+ fik In 


\ — El 
p. 


I Plp2 
P 2 „2 ’ 


Pt 


( 2 ) 


with A and fi being the Lame constants, ps the saturated dislocation density, 
and k and y the constants of the material. The first two terms in ([^ represent 
the contribution to the energy of crystal by the elastic strain. The third 
and fourth terms are the contributions to the energy from the two separate 
groups of excess dislocations (Berdichevsky, 2006b). Finally, the last term 
takes the interaction of dislocations between these two groups into account. 
For small up to moderate dislocation densities the logarithmic terms can be 
approximated by the Taylor series as follows 

In-- ^ — + -P-. 

^ Ps Ps ^ Ps 

We shall use further only this approximation. 

With ([^ and ([^ the total energy functional of the bent beam becomes 


i‘h/2 


I = a 


1 1 1 r 

2 + Uy^yf + 2(y9i + -(32 Sm 2(^2)' 


'-h/2 L 

W . W . .. .2 1 


+ p{uy^y - -f3i sin 2(pi - -^2 sin 2(^2) + ^p{ux,y + - Pi cos 2ipi 


P 2 COS 2p>2Y + pk 


cos</9i + Pi^y sin (fil 


bps 


, 1 , {Pl,xCOSipi + Pi^ySiULfiY 

2^ yv. 


pk 


\P2,x COS<y92 +/?2,ySinV32| 


bps 


, 1 , (/52,xCOSV 32 +/^2,j/Sm(^2) 

+ - W~2 - 

2 h^pl 


+ kpx 


\Pi^x cos(pi + Pi^y sin(pi||/32,x C 0 S(P 2 + P 2 ,y sin(p 2 | 


h^pl 


dxdy 


i-h/2 


'-h/2 


ryUx\x=Ldy. (3) 


The last term in (|^ is the work done by the linearly distributed traction —ry 
acting at the boundary x = L. If the dissipation caused by the dislocation 
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motion is negligible, then the true displacements Ux, Uy and plastic slips /3i, 
f52 in the final state of bending minimize energy functional (|^ among all 
admissible displacements and plastic slips satisfying the kinematic boundary 
condition. By the hnal state of bending we means the equilibrium state which 
is established at hxed loading condition after dislocation nucleation and after 
the movement of dislocations toward their equilibrium positions is hnished. 
Thus, the whole process of dislocation nucleation and dislocation dynamics 
cannot be described by this theory. The bending moment M = rah^/12 
is regarded as a control parameter, so one can study the evolution of the 
dislocation network in terms of M. 

If the dissipation due to dislocation motion cannot be neglected, energy 
minimization should be replaced by the variational equation (Sedov, 1965) 

pL 


Sl + a / {^6f3i +^6f32)dxdy = 0. 

Jo J-h /2 opi dp 2 


(4) 


The last integral in this equation describes the energy dissipation due to 
dislocation motion, where the dissipation function D{fii,fi 2 ) is assumed to 
depend only on the plastic slip rates. We shall consider the simplest rate- 
independent theory for which 

-D(A,/^2) = K{\iJl\ + |/^2|), 

with K being the critical resolved shear stress. Then, provided the signs of 
and ^2 do not change during the evolution of fii and /92, the variational 
equation (|^ reduces to minimizing the following “relaxed energy” functional 

pL ph/2 r ^ ^ 

dd = a / -X{Ux,x + Uy,yY + p{ux,x +7:fiism2ipi +-/32sin2ip2fi 

Jo J-h /2 I I 

-h p{uy^y - sin 2ipx - ^[^2 sin 2 ip 2 )^ + ^p{ux,y + Uy^x - fii cos 2ipi 


- [32 cos 2 ip 2 Y 


\l3i^x cos<y9i -F/3i,ySin<y9i| 


hps 


, 1 ^ {[Jl^xCOSPlX- fil^ySinpiY 

2^ m 


fik 


1/52,X COS(p 2 + [32,y sin ip 2 


bps 


1 , {(32,X cos p >2 + (32,y sin <^ 2 )^ 


kpx 


\(3i^x cos Pi + (3i^y sin pi \\(32,x cos p 2 + (32,y sin p 2 \ 


b‘^p\ 


K sign{(3i)(3i K sign(/? 2)/?2 


ph/2 

dxdy + a / TyUx\x=Ldy. (5) 

J-h/2 
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So, the true displacement and plastic slip fields in the final equilibrium state 
of plastic bending minimize the “relaxed” energy functional Id among all 
admissible displacements and plastic slips. Finally, if /3i = 0 and /32 = 0, 
then the plastic slips are frozen, while the displacements should be found by 
minimizing ([^ with the frozen plastic slips. Let us mention here also the 
more complicate model taking into account the interaction of dislocation in 
dissipation leading to the latent hardening. In such model we assume the 
rate-idependent dissipation function in the form 

$ 2 ) = K{\$i \ + \fi2 \ + 

In this paper we shall not dwell on this model. 


3. Energy minimization 

We first focus on the minimization problem (|^. It is convenient to intro¬ 
duce the following dimensionless variables and quantities 


X = bpsX, y = bpsy, h = bpsh, L = bpsL, r = 

VjX bpstlx-, Uy bpgXiy^ 1 , y , 

pa p 


T 


pbps' 


( 6 ) 


which simplify the minimization problem. Now the energy functional can be 
rewritten in the dimensionless form as follows 


/ = 


lo J-h/2 _ 


^'y{ux,x + Uy^yY {ux,x + sin 2lpi + sin 2 Lp 2 f 


+ sin 2ipi - sin 2 ip 2 f {ux,y + Uy^x - fix cos 2ipi 

- fi2 COS2(^2)^ + k\l3i^xCOS(pi + fii^ySmifil + ^k{/3i^x cos(pi fii^ySimpiY 
+ k\fi 2 ,x cos lf 2 + /32,y sin If 2 \ + ^k{/32,xCOSip2 + /32,ySinf2Y + kx\fii,x cos If I 

j‘h/2 

dxdy+ / TyUx\x=Ldy. (7) 
J-hl2 


+ fii^y sin fi 11/32,X cos f 2 + fi 2 ,y sin f 2 \ 


The bar over the dimensionless quantities is dropped for short because in the 
subsequent analysis we shall deal only with them. 
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Energy functional Q contains a small parameter h, so it can be re¬ 
duced to 1-D energy functional by the variational asymptotic method (see 
(Berdichevsky, [1983 ; |1999 )). For this purpose let us introduce the rescaled 
coordinate ( = y/h, ( E (—1/2,1/2) to £x the domain over the thickness in 
the passage to the limit h —)■ 0. With this new variable ( the small parameter 
h enters the functional explicitly through the formulas 

Since the boundary condition aX x = L does not influence the inner asymp¬ 
totic distributions of the displacements and the plastic slips over the thick¬ 
ness, we put r = 0 for a while. At the hrst step of the variational-asymptotic 
procedure we keep the asymptotic principal terms in ([^ to obtain 

k k k 

+ ^(A,csinv^i)^ + -\/32,(Sm(p2\ + — (/? 2 ,c sin 


+ -j^\/3i,(sm^i\\/32,t; smip2\ 


dxd(. 


Functional Jq is positive dehnite, so its minimum must be zero and is achieved 
at = 0. For the bending states which we are 

interested in let us set at this step = 0 and /3i = /da = 0 to obtain 

Ma: = 0, Uy=v{x), /?1 = ,da = 0. 

At the second step, we fix v{x) and seek the minimizer in the form 

M^ = <(a;,C), Uy = v{x)+u'y{xX), /^2 =C), 

with u'^{xX),Uy{xX)y P 2 i^X) being the correction terms which are 

small compared to v{x). Substituting these formulas into ([^ and then 
keeping the asymptotically leading terms containing m/(x,/'), Uy{xX), and 
(3'^{x, C), we obtain 

r r 


" 1/2 [1 1 11 U 

+ 7:2 «c)^ + b(F<c + 


' 1/2 


h2 


k 


k 


2'-h 
k 


h' 

Xk, 


+ ^(/^pcsin/3i) + j^\l32^^smip2\ +^{/32,^smip2) + 


2h2 


h2 


1/^2,/ siiKpal 


dxd( 









Since functional Ii is also positive definite, its minimum is again zero and is 
achieved at 


= 0’ “kc = = 0, /?2,c = 0- 

This implies 

u'y = 0 , = -hCv,^, = 0 , (32 = 0 . 

At the third step of the variational-asymptotic procedure, we look for the 
minimizer in the form 

Ux = -hCv^x + u^{x, c), Uy = v(x) + Uy{x, C), (8) 

/3i = /?"(a;,C), /32=/3'^{x,C), 

with u'(.{x, (), Uy(x, 0, (3'({x, (), and (32{x, () being the correction terms which 
are assumed small compared with those found at the second step. By redehn- 
ing v(x) if required, we may put the following constraints on these correction 
terms 


(“P = (^ 2 ) = 0) where (•) 



(9) 


Constraint (Uy) = 0 means that v(x) is the average deflection (over the thick¬ 
ness) of the beam. Substituting ([^ into ([^, then keeping the asymptotically 
leading terms containing u'(.(x,(), Uy(x,(), /3'((x,(), and /32(x,() to get 


l9 = h 


l>L i>l/2 

'0 Jl/2 . 
1 


1 1 (3” (3” 

-7(-hCn,,,r+^M"c)^+(-^C^^,xx+^sin2(pi-F^sin 2<^2) 


1.1 


+ sin2(pi - -(32 sin 2 (^ 2 ) + “ /^i cos2(pi - /?2 cos 2 (^ 2 ) 


2"h 


k 


+ ^l/^iysinv^il + — sin v?i) + -|/32,^sinv72| + —{^ 2 ,<:^^^^ 2 ) 


2h2 


h' 


I/O// • 11 /P" • I 

^|Pi,c sm(pi||/t2,c sm(p2| 


2h2 

Jff 


dxd(. 


Functional I 2 can be reduced to a functional depending only on jd” and /d!^. 
Indeed, hxing hrst f3'{ and (d'^ and varying this functional with respect to u'(. 
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and u'' then nsing the natnral bonndary conditions at C = ± 1 / 2 , we obtain, 


-(7 ± 2 )m"^ = ± P'l sin 2 v 3 i ± sin 2932 , 


-m" ,. = f3” cos 2931 ± (32 cos 2 ip 2 - 


( 10 ) 


After hnding n" and u” according to these eqnations, we snbstitnte (10) into 
0 and change ( back to y. Since the fnnctional does not contain and 
(32^, the thickness problem rednces to minimizing the following fnnctional 
with respect to /3”{x,y) and /32{x,y) 

ph/2 


h = 


’-h/2 . 


K 


-{-2v^xxy ± /3/sin2(pi ± sin 2932 ) ± k\l^'ly sin(pi| 


+ ^(/^i,j/Sin(pi)2 + k\PlySixnp2\ + ^{P'lySixnp2f 


+ xk\l^i,y sin</5i||/52,j^ sin 932 1 


dy, 


where k = 


1 


2(l-u)- 


Solntion of this minimization problem enables one to 


hnd 13” {x, y) and /?2 {^1 v) reqnired for the rednction to 1-D theory. The 
analytical solntion of the thickness problem can be fonnd in the special case 
of symmetric slip systems with the angles ipi = ip & ( 0 , 7 r/ 2 ), 752 = tt — </ 3 . 
Obvionsly, f3” = f3 = —^2 in this case dne to the symmetry of the problem. 
Thns, fnnctional Ja simplihes to 

I>hl2 

h= / 




>-h/2 


^i-‘2v^xxy ± 2/3sin 293 )^ ± 2fc|/3,j^ sin<y3| ± k{l ± x)(/3,y sin93)^] dy. 

_ ( 11 ) 

Fnnctional (pTj) is similar to the fnnctional obtained in (Le and Nguyen| 


2013), so we present the solntion of the thickness problem withont detailed 


derivation. The minimizer /3{x,y) is 

pm (a:) 

P{x,y) = < /3o{x,y) 

[-P{x,-y) 

Fnnction /3o{x,y) is fonnd to be 

1 


for y e {-h/2,-l{x)/2), 
for y e (-/(a;)/2,0), 
for y e (0, h/2). 


( 12 ) 


Mx,y) = - 


',xx 


sin 2(p rj 


W 


cosh ^ 


sinh rjy 
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where 


r; = 2. 


2k 


k{l + x) 


cos if. 


Function (3m{x) is given by 


+ tanh ^ 
sin Iiprj \ 2 2 

Finally, the length l{x) should be found by solving the following transcen¬ 
dental equation 

—^-1-tanh ^)(h—/)] sin2(^ —fcsign(/3o,y) sin(^ = 0. (13) 

If the length I is small, then the sign of v^xx coincides with the sign of 
evaluated right from the point y = —1/2. For dehniteness let (3o^y{—l/2 + Q) > 
0 so that v^xx > 0. Note that, if the curvature of the beam is constant, then 
l{x) does not depend on x and remains constant over the whole length of 
the beam. In this case (13) can be regarded as the equation for v^xx once I 


is known. By integrating (10) in our special case and taking the constraints 
into account, we obtain the displacements Ux and Uy in the form 


'^x '^,xyi 

= + + (3{x,0d^- ^^sm2ip, 

wheTeA = {J^/3{x,0d0- 

Having found the solution of the thickness problem, let us now substitute 
the above formulas for the displacements together with the expression of f3 
into the energy functional ([^. Keeping the asymptotically principal terms 
and integrating over the thickness we obtain the one-dimensional functional 

J[v{x)]= [ ^{v^xx) dx - Mv^x\x=l, (14) 

Jo 

where the bending energy density reads 

<h(a;) = cica^ -I- C2UJ. 
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In these formulas we use uj to denote the curvature v^^x, M = rh?/VI is the 
resultant moment, and 

/•-h 2 /-o 


Cl = 2 
C 2 = 2 

with 


/ (2go-22/)^d2/ + 2 

l-h/2 ^ 

f = 

/_//2 sm293 COS</9 


- 1/2 


-{2q-2y) +-- —{q,y) 

2 4 cos"' ip 


2k 


?o, 


(l{y) = - 

V 


vy - 


cosh ^ 


sinh rjy 


dy, 

(15) 

(16) 


and go = Note that, as the coefficients Ci and C 2 depends on the 

curvature u through I, the energy density is not quadratic with respect to u. 
Besides, if q{y) = 0, then ci = nh^/G and C 2 = 0, so the obtained functional 


reduces to the classical 1-D functional of the elastic beam as expected (Le 


1999). 


Varying functional (|14|) with respect to the deflection v, we obtain the 


differential equation of bending 


'^,XX O 5 

subject to the boundary conditions 


m = 


du 


(17) 


(18) 


(19) 


11 ( 0 ) = 0 , vJO) = 0 , 

m{L) — M = 0, m^x{L) = 0. 

Equation 0. together with the conditions ([ 1 ^ 2 , implies that 

/ \ ^ dci 9 dco 

Tfiioj^ = 2c\ijj V ——ijj C 2 ~—00 = M. 

du doj 

Since the bending moment m is independent of x, the curvature must also be 
constant over the length of the beam. Consequently, I and /Sq are independent 
of X and the upper and lower layers of the beam are dislocation-free. Together 
with (13), this equation determines the moment-curvature curve during the 
plastic bending. To plot this curve let us compute the derivatives of ci and 
C 2 with respect to uj 

dci dci dl dc2 dc2 dl 
duj dl du^ du dl du 
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From (13) we find that 


duj 4:kri tanh {rj {h — 1) tanh ^ + 2) 


dl 


K COS (p{h — ly {r]{h — /) + 4 tanh 


2 • 


This formula, together with 10). enables one to determine dci/doj and 
dc 2 /d(jj^ required for plotting the moment-curvature curve. 

The threshold value of curvature at which dislocations begin to nucleate 
is calculated by letting I go to zero in (13). This yields 


2k 




cos if 


The threshold value of moment can be computed from ( |19[ ). Taking into 
account that, at / = 0, ci = C 2 = 0, while dci/doj = dc 2 /du = 0, we 

get 


= 


2kh 


3 cos (f 

Returning to the original physical quantities according to ([^ we have 

2fikah 


M™, = 


Tp. 


12 


3bps cos ip ’ 


which exhibits explicitly the size effect. Besides, M^n goes to inhnity as 
the shear modulus p goes to inhnity, so this coincides well with the result 
obtained by Ashby (1970). Thus, if M < Men, then /3 = 0, so no dislocation 
is nucleated and we have purely elastic solution. The plastic yielding begins 
at M = Men- 

Combining the elastic and plastic responses together, we get the moment- 
curvature relation in the following form 


M = 


/ 


Kh^ 


bJ 


for M < Me- 


\2cia; + + C 2 -1- for M > Men- 

Knowing the curvature u: from M, we integrate the equation v^xx = ^ and 
use the boundary conditions (18)i to obtain the dehection of the beam 


vix) = -ux^. 
^ ’ 2 
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Figure 2: Function Men{v) 
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Figure 3: Function fi{y) for M = 0.0001 and (a) (p = tt/S, (b) p = tt/O 


4. Numerical simulations 


In order to simulate the minimizer numerically, we choose h = 0.1, L = 
10, u = 0.25, k = 10“"^, and y = 0.8. Fig. @ shows the plot of energetic 
threshold of the bending moment Men as function of the angle ip. It can 
be seen that the threshold moment becomes infinite as ip goes to 7r/2 and 
has a local minimum at (p = 0. It must be emphasized that this result is 
obtained under the assumption that the crystal admits two active symmetric 
slip systems. So, if (p = 7r/2 the slip systems merge into one that remains 
inactive during the whole process of bending because of its vanishing resolved 


shear stress (cf. with (Le and Nguyen, 2013)). However, if the slip systems 


are oriented unsymmetrically, then the one with largest Schmid factor will 
be activated at a certain hnite bending moment. 
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For M > Men the plastic slip becomes non-zero. Fig. |^show the plots of 
plastic slip /3{y) for M = 0.0001 at two different orientations of slip systems. 
The plastic strain vanishes on the middle line of the beam as expected and 
reaches its maximum and minimum at the traction-free faces. There are 
two small boundary layers near these faces where the plastic slip remains 
constant. 


P 



Figure 4: Dislocation density for M = 0.0001 and (a) ip = tt/S {I « 0.097), (b) ip = 7r/6 
[I « 0.099) 


V 



Figure 5: Deflection of the beam 

On Fig. where the dimensionless dislocation density p{y) = sin ip is 
plotted for M = 0.0001 and for two different angles, it is seen that the excess 
edge dislocations of the same sign of each group are concentrated in the 
middle of the beam thickness, with maximum dislocation density achieved 
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at 1 / = 0. Although no obstacle exists on the middle line, the repulsive forces 
between dislocations of the same sign prevent them from colliding. Thus, the 
high concentration of dislocations can be regarded as the dislocation pile-up 
against the middle line. The dislocation free zones are y G (—h/2, —1/2) and 
yeil/2,h/2). 

The dimensionless deflection of the beam, v{x), is shown in Fig. for 
M = 0.0001 and (a) ip = tt/S, (b) ip = tt/O. Since the curvature of the beam 
is constant, the thickness of the dislocation zone I does not depend on x. 



Figure 6: Moment-curvature curve for: (a) ip = 7r/3, (b) p = 7r/6 

On Fig. 1^ we show the moment-curvature curve for <^9 = vr/S and ip = vr/O. 
Up to the threshold moment M^n (corresponding to point A on this figure) 
the moment-curvature curve is a straight line corresponding to the linear 
elastic beam theory. Then the curve becomes non-linear and increasing as M 
increases and I increases from zero to h/2. This nonlinear portion describes 
the work hardening due to the dislocation pile-up against the middle line of 
the beam. If the bending moment is increased from zero up to the moment 
Mb corresponding to B (the loading case) we follow the moment-curvature 
curve from O through A to B. Now, if we unload the beam by decreasing 
the bending moment from Mb to zero, the curvature and moment follow the 
same path BAO (in the inverse direction), and at the end of the unloading 
process no residual curvature of the beam is observed. 
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5. Non-zero dissipation 

If the dissipation cannot be neglected, the problems reduce to minimizing 
the relaxed energy functionals ([^. The signs of the dissipation terms in this 
functional depend on whether /3„ > 0 or < 0, a = 1, 2. However, for the 
beam bending it is easy to see that both cases occur simultaneously during 
the plastic deformations. Indeed, from the elementary beam theory (and 
also from the previous simulations) we know that, as the bending moment is 
increased (loading), /3i > 0 for y > 0 and /9i < 0 for y < 0, while /92 < 0 for 
y > 0 and /32 > 0 for y < 0. In contrary, if the bending moment is decreased 
(unloading or loading in the opposite direction), f3i is either frozen or /3i > 0 
for 1 / < 0 and /di < 0 for ?/ > 0 while the opposite happens to /32. Since 
these functionals differs from each other in the loading and unloading case, 
the case study must be done separately. 

The case study becomes simpler for the symmetric slip systems with 
(pi = (p G (0,7r/2), ip 2 = 71 — (f. Due to the symmetry of the problem, 
we have (3i = {3 = —/? 2 - Therefore ii'sign(/di)/3i = ii'sign(/d 2 )/? 2 - According 
to the above property, both terms can be replaced by Ksignyfi. We use 
again the dimensionless quantities (|^ and write the functional ([^ in this 
special case in the form 


Id = 


-jl {Ux,x + Uy^yf + (Ua:,* + (3 siu 2Lpf + {Uy^y - 13 sin 2Lpf 

■^2 ^ 2 
+ - (Mx,?; + My,x) + k\l3 ,x COS (f + ^^y siu (y9| + -k {(3^^ COS If + (3^y sin ip) 

1 2 
+ k\(3^^ cos ip - [3^y sin (p| + -k {/3^x cos p) - I3^y sin ip) 



+ kx\/3,x cos ip + (3,y sin v?! ■ |/3,a, cos p - (3^y sin v?] + 2e sign?/ j3 

rh/2 

+ / Tyux\x = Ldy, 

J-h/2 


dxdy 


( 20 ) 


where e = K/y. As compared to the previous functional ([^ the additional 
term 2e sign y (3 does not belong to the asymptotically principal terms. There¬ 
fore, up to the second step of the variational-asymptotic procedure this term 
does not have any influence on the inner asymptotic expansion. At the third 
step, we are looking for the minimize in the form ([^ such that the constraints 
(|^ are fulhlled. Fixing (3” and minimizing the relaxed energy with respect 
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to m" and Uy we find them in the form (10). Then, the fnnctional rednces to 


r‘h/2 r 


h = 


^ {-2v^xxy + 2f3"sm2(pf + 2k\(3''ySm(p\ + k{/3”y sin (pf 


J -h12 

+ xk {(3”y sin ipY + 2e sign?/ /3" 


dy 


Now the term 2esign?/ /3" shonld be kept because it has the same order as the 
cross term —4:Kv^xxyld"&in2ip. Up to an unessential constant we may rewrite 
this functional as 


h = 


rh/2 

l-h/2 


'x f -‘2v,xxy + + 2/5" sin 2(p ] + 2fc|/5" sin (p| 


+ + xW'ySimpf 


K sin 2(p 
dy. 


( 21 ) 


Functional (21) is similar to the functional obtained in (Le and Nguyen| 


2013) for the case of non-zero dissipation, so we again present the solution of 
the thickness problem without detailed derivation. We use the same Ansatz 


(12) for /3 in this case. Then for function jdoiy) we hnd 


^ 0 ( 2 /) = 


1 

sin2(y9 r] 
e 


7]y 


cosh ^ 


sinh rfy 


2fi;sin^ 2p> 
Computing /5o(2/) ai y = —1/2 we get 


(1 — cosh rjy — tanh ^ sinh r/?/) . 


( 22 ) 


/J™ = (-'2+ tanh ^2) + y 

sm 2<pri \ 2 2 J 2k sin^ 2ip 

The transcendental equation for I reads 

f^v^xx[^ih^ - Y + tanh ^)(h - 

- I sign(/3i,j,) sin (p = 



cosh ^ 


- /)] sin2(p 

0 . 
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Knowing jS, we can now determine the displacement field in accordance 
with (10). Then, substituting this field together with (]{y) from (12) and (22) 
into the energy functional (20), keeping the asymptotically principal terms 
and integrating over the thickness we obtain the 1-D functional 


Jd[v{x)]= ^{v^^^)dx - 

Jo 

with the bending energy density 

$(0;) = cio;^ + C2UJ. 
The coefficient ci and C 2 are given by 


(23) 


Cl = 2 

C2 = 2 

-2y)) 


r-l/^ rO 

I (2qo-2yfdy + 2 I 

-h/2 ^ 

0 r ^ e ^ ^ 


J-112 
,/c(l + x) 


^ in 0 ^ 2 , k{'^ + x), .2 , 


-1/2 [COSV? 
dy + 2 I 


K sin 2(p ^ 4 cos^ cp 


l>-i/2 

-h/2 L 


sin2(p 


(po - l)(2go - 2p) 


dy, 


with 


1 smhpp, 

Qijj) = -ijiy - ^2 

h cosh f 


Til 

p{y) = 1 — cosh rjy — tanh — sinh pp. 


andpo = p(~V2), do = <?(—V2)- Note that the coefficients Ci and C 2 depends 
on the curvature cn, so the energy density is not quadratic with respect to oj. 


Varying functional (23) with respect to the deflection v{x) we get the same 
equilibrium equation and boundary conditions as in Section 4. Consequently, 


the moment-curvature relation (19) remains valid if ci and C 2 from above are 


substituted. It turns out that the straight line M = nh^u/S corresponding 
to the purely elastic solution does not intersect the moment-curvature curve 
corresponding to the solution with dislocations at / = 0 for the case with 
non-zero dissipation. The threshold value for the curvature must therefore 
be calculated by solving the system of equations 


= uJd, rn(ld) = Kh^Ud/3, 


(24) 
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where and ujd are regarded as unknowns]^ Using (19) we determine also 
the dissipative threshold value for the bending moment Computing the 
curvature from the given moment, we hnd the deflection by integrating the 
equation 

Consider now the following loading/unloading process. We hrst bend the 
beam slowly and successively by increasing the moment until some maximum 
value M* > Md such that at the end of the loading process the plastic slip 
becomes /3*. Then we unload the beam by reducing the bending moment 
back to zero. Assuming that dislocations are frozen during this unloading 
process with (3 = (5^ and /d = 0 , we have to find the displacements of the 
beam by minimizing the energy functional 


Ia = 


rh/2 

'-h/2 L 


{Ux,X + Uy^yf + A (3 SlU 2^pf + {Uy^y - /3 SlU 2 " 


1 2 1 
+ 2 + ^y,^) + ^1 Aa: cos (p + f3^y sin (p| + -k {j3^^ cos (p + j3^y sin (p) 

1 2 

+ cos (p — f3^y sin (p| + -k {j3^x cos (p — (3^y sin ip) + kx\(3,x cos ip 


+ f3^y sin 1 ■ \f3^^ cos ip - (3^y sin 9 ?| + 2e signy (3 


rh/2 

dxdy+ / TyUa;\x=Ldy. 

J-h/2 

(25) 


where, now, f3^ is hxed and does not subject to the variation. Since functional 
(25) is functional ([^ with cpi = n — (p 2 = (p and /9i = —(32 = (3^, exactly the 
same variational-asymptotic analysis as in Sections 3 can be repeated leading 
to the following displacement held 




u 


y = v{x) + ^^{y -i2)Uxx + ;^ 


(3^ (x, ^)d^ - ) sin 2ip, 


where A* = (3^{x,^)d^). We substitute this displacement held into func¬ 

tional (25) and keep the asymptotically principal terms. Integrating over the 


thickness and neglecting the terms containing the known function (3^ only. 


^Formula for uJd and Figure 7 plotting in (Le and Nguyen 20131 are not correct 
and should be amended in accordance with (24). 
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we get 1-D functional (23) with different coefficients Ci and C 2 


l>o r-lfi 

ci = Kh^/6, 02 = 2 1 Ku^q^{y){-4y)dy+ 2 Ku^qo^{y){-4y)dy, 

J-l/2 J-h/2 

where ca* is the curvature corresponding to M* and 

(l*iy) = - Ivy -^^^sinh?]?/) , go* = g*(-/*/2), 

V \ cosh^ J 

U being the length corresponding to (5^. Thus, during this unloading process 
the moment-curvature relation takes the form 

m = 2ciUJ -f C 2 = M, 

so that, at the end of the unloading when M = 0, the residual curvature is 

C2 

Ui = -. 

Cl 

Note that, although the bending moment after unloading is zero, the elastic 
strain is not, so there is still some eigenstress in the beam. 


6. Numerical simulation 


In order to simulate the minimizer numerically, we choose the same nu¬ 
merical values for h, L, u, k, x Section 4. For e = K/fi we choose the 
value e = 0.0001. In order to hnd the dissipative threshold for the bending 
moment we need to solve equations (24) with respect to Ud and Id- Then 


Md = nh^Ud/S. Numerical simulations give for instance 

ujd = 0.0413739, Id = 0.0192525 ^ Md = 0.00001488 for (p = 7r/3, 

and 


Ud = 0.0669535, Id = 0.0198977 ^ Md = 0.91942 x 10“® for cp = 7r/6. 

Thus, if we increase M steadily, then at the onset of plastic yielding the 
plastic slip and the total number of dislocations turn out to be hnite. This 
could be explained physically in the following way: in the presence of non¬ 
zero dissipation dislocations should accumulate enough to give the crystal 
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Figure 7: Function j3{y) for M = 0.0001 and (a) tp = tt/S, (b) ip = tt/G 


lower “relaxed” energy as compared to the energy of pnrely elastic bending. 
This means also that the states with smaller number of dislocations present 
energy barriers to be overcome, and to have the smaller “relaxed” energy the 
crystal should jump to the certain hnite dislocation density. 

For M > Md the plastic slip becomes non-zero. Fig. shows the plots of 
plastic slip /3{y) for M = 0.0001 at two different orientations of slip systems. 
The plastic strain is zero on the middle line of the beam and reaches its 
maximum and minimum at the free faces. There are two small boundary 
layers near these faces where the plastic slip remains constant. 



Figure 8: Dislocation density for M = 0.0001 and (a) p = tt/S (Z « 0.097), (b) p = 7r/6 
(Z « 0.099) 

On Fig. 1^ where the dimensionless dislocation density p{y) = f3^y sin (p is 
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plotted for M = 0.0001 and for two different angles, it is seen that the excess 
edge dislocations of the same sign of each group are concentrated in the 
middle of the beam thickness, with maximum dislocation density achieved 
at ?/ = 0. Although no obstacle exists on the middle line, the repulsive forces 
between dislocations of the same sign prevent them from colliding. Thus, the 
high concentration of dislocations can be regarded as the dislocation pile-up 
against the middle line. The dislocation free zones are y G {—h/2, —1/2) and 
ye{l/2,h/2). 


V 



Figure 9: Deflection of the beam 


The deflection of the beam v{x) is shown in Fig. [^for M = 0.0001 and 
(a) ip = vr/S, (b) ip = vr/O. Since the curvature of the beam is constant, the 
thickness of the dislocation zone I does not depend on x. 


On Fig. 10 we show the moment-curvature curve for (p = tt/S and (p = 
vr/O, respectively. Up to the threshold moment Md (corresponding to point 
A on this hgure) the moment-curvature curve is a straight line corresponding 
to the linear elastic beam theory. Then the curve becomes non-linear and, 
in contrast to the zero dissipation case, also non-monotone as uj increases 
as can be seen in curve (b). At the onset of the plastic yielding we can ob¬ 
serve a small softening behavior of the moment-curvature curve. The further 
monotonously increasing portion of the moment-curvature curve describes 
the work hardening due to the dislocation pile-up against the middle line of 
the beam. If the bending moment is increased from zero up to the moment 
Mb corresponding to B (the loading case) we follow the moment-curvature 
curve from O through A to B (with a possible jump in uj to pass the soft¬ 
ening portion). If the beam is then unloaded, the moment-curvature curve 
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Figure 10: Moment-curvature curve during loading and unloading for: &) ip = tt/S, b) 
ip = 7r/6. 


becomes a straight line BC with the same slope like that of OA as shown in 


Fig. 10 


7. Discussion and outlook 

In this paper the one-dimensional theory of bending of a single crystal 
beam with two active slip systems has been developed. The one-dimensional 
energy of bending can be found after the solution of the thickness problem. 
The latter has been solved for the symmetric slip systems. The threshold 
bending moment exhibiting the size effect has been found for the case without 
and with dissipation. We have found also the dislocation density and the 
moment-curvature curves at loading and unloading. It would be interesting 
to verify the obtained results with the measurement of excess dislocation 
density in a bent beam by the method of orientation imaging microscopy 
developed recently in ( |Kysar et ah 2010| . The specimen could be for instance 
the fee single crystal beam having the axis parallel to [ilO]-direction and bent 
about the [110] direction. For the unsymmetrically oriented slip systems the 
thickness problem can only be solved numerically. This will be addressed in 
our forthcomming paper. 
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